# Economic Risk within the Household and Voting for the Radical Right
# Tarik Abou-Chadi and Thomas Kurer
# World Politics

## REPLICATION FILES

# Descriptives (Figures 2 and 3)

# Preparation
rm(list = ls()) 
dev.off()
cat("\014")  
# globals
options(scipen=999)

# packages
if (!require("pacman")) install.packages("pacman") 

pacman::p_load(tidyverse, broom, hrbrthemes, plm, estimatr, sandwich, lmtest, AER, lfe, huxtable, margins, readstata13, texreg, reshape2, readxl, xtable, haven, ggrepel, magrittr)

# ggplot theme
theme_set(theme_bw() + theme(text = element_text(size=16)))

# set paths (needs to be adjusted before running script)
path <- "..."

# load data

df <- read.dta13(paste0(path, "riskradright.dta"))

# unemp risk by occupational group

df$isco_1d <- ifelse(df$year<2012, df$isco88_1d, df$isco08_1d)
df$isco_1d[df$isco_1d==0] <- NA

df$isco_1d_label <- as.factor(df$isco_1d)
levels(df$isco_1d_label) <- c("Managers","Professionals","Technicians", "Clerks", "Service/Sales", "Skilled Agricultural", "Craft", "Machine Operators", "Elementary Occupations")

df$female_label <- as.factor(df$female)
levels(df$female_label) <- c("male", "female")

df %>% filter(!is.na(female_label) & !is.na(isco_1d_label)) %>% group_by(isco_1d_label, female_label) %>%
  summarise(unemprisk=mean(unemp_isco_1d, na.rm=T)) %>%
  ggplot(aes(x=isco_1d_label, y=unemprisk, group=female_label, fill=female_label)) + geom_bar(pos="dodge", stat="identity") +
  xlab("") + ylab("Occupational Unemployment Risk") + scale_fill_grey(start=.4, end=.8) +
  theme(text = element_text(size=16), axis.text.x=element_text(angle=45, hjust=1, vjust=1, size=14), legend.title=element_blank(), legend.justification = c(0, 1), legend.position = c(0, 1), legend.box.margin=margin(c(5,5,5,5)))

ggsave(paste0(path, "unemprisk_isco_female_compositional.pdf"), height=6, width=6)


# household risk constellations

heat <- df %>% dplyr::select(cntry, year, gndr, unemp_isco_2d, punemp_isco_2d) %>% dplyr::filter(!is.na(unemp_isco_2d)) %>% dplyr::filter(!is.na(punemp_isco_2d))

heat %<>% group_by(cntry, year) %>% mutate(perc_ur=ntile(unemp_isco_2d, 10), perc_pur=ntile(punemp_isco_2d, 10)) %>% ungroup()

heat$combo <- paste(heat$perc_ur, heat$perc_pur, sep="-")

heat_all <- heat %>% mutate(one=1, total=sum(one)) %>% group_by(combo) %>% mutate(freq=sum(one)) %>% dplyr::filter(row_number(cntry) == 1) %>% ungroup() %>% mutate(perc_freq=ntile(freq, 100))

heat_male <- heat %>% filter(gndr==1) %>% mutate(one=1, total=sum(one)) %>% group_by(combo) %>% mutate(freq=sum(one)) %>% dplyr::filter(row_number(cntry) == 1) %>% ungroup() %>% mutate(perc_freq=ntile(freq, 100))
heat_female <- heat %>% filter(gndr==2) %>% mutate(one=1, total=sum(one)) %>% group_by(combo) %>% mutate(freq=sum(one)) %>% dplyr::filter(row_number(cntry) == 1) %>% ungroup() %>% mutate(perc_freq=ntile(freq, 100))


ggplot(data = heat_all, mapping = aes(x = perc_ur, y = perc_pur, fill = perc_freq)) +
  geom_tile() +
  scale_fill_gradient2(low = 'white', mid = 'gray90', high = 'black') + 
  scale_x_continuous(name="Unemployment Risk Respondent (Decile)", breaks=seq(1,10, 1)) +
  scale_y_continuous(name="Unemployment Risk Partner (Decile)", breaks=seq(1,10, 1)) +
  ggtitle("All Respondents") +
  theme(legend.position="none", text = element_text(size=16))
ggsave(paste0(path, "heat_all.pdf"), height=6, width=6)
  

ggplot(data = heat_male, mapping = aes(x = perc_ur, y = perc_pur, fill = perc_freq)) +
  geom_tile() +
  scale_fill_gradient2(low = 'white', mid = 'gray90', high = 'black') + 
  scale_x_continuous(name="Unemployment Risk Respondent (Decile)", breaks=seq(1,10, 1)) +
  scale_y_continuous(name="Unemployment Risk Partner (Decile)", breaks=seq(1,10, 1)) +
  ggtitle("Male Respondents") +
  theme(legend.position="none", text = element_text(size=16))
ggsave(paste0(path, "heat_male.pdf"), height=6, width=6)
  
ggplot(data = heat_female, mapping = aes(x = perc_ur, y = perc_pur, fill = perc_freq)) +
  geom_tile() +
  scale_fill_gradient2(low = 'white', mid = 'gray90', high = 'black') + 
  scale_x_continuous(name="Unemployment Risk Respondent (Decile)", breaks=seq(1,10, 1)) +
  scale_y_continuous(name="Unemployment Risk Partner (Decile)", breaks=seq(1,10, 1)) +
  ggtitle("Female Respondents") +
  theme(legend.position="none", text = element_text(size=16))
ggsave(paste0(path, "heat_female.pdf"), height=6, width=6)

# how many are in those cells?

heat$extreme <- ifelse(heat$combo%in%c("1-10", "2-10", "3-10", "1-9", "2-9", "3-9", "1-8", "2-8", "3-8", 
                                       "10-1", "10-2", "10-3", "9-1", "9-2", "9-3", "8-1", "8-2", "8-3"),
                                       1, 0)

table(heat$combo, heat$extreme)

sum(heat$extreme)/nrow(heat)
# 12%

offdiagonal <- heat %>% filter(cntry!="LU" & cntry!="BE" & cntry!="PT" & cntry!="IT" & cntry!="IE" & cntry!="ES") %>% mutate(one=1)  %>% group_by(cntry) %>% summarise(total=sum(one), extreme=sum(extreme), share=(extreme/total)*100)
print(xtable(offdiagonal, caption="Off-Diagonal Population Share", label="fig:offdiag", digits=c(0,0,0,0,1)), file = paste0(path, "offdiagonal_share.tex"))

